Este R Markdown recoge el enunciado de la práctica de la asignatura de redes sociales.
El objetivo es analizar un grafo, que se provee como fichero en el mismo paquete que este enunciado. En este fichero, encontramos solamente dos columnas, correspondiente a una interacción entre dos nodos de la red. Esta red está formada por distintos individuos que tienen contactos cara a cara durante un período de tiempo.
A continuación, dividimos la práctica en apartados, con una breve descripción de qué debe contener cada chunk de código donde el alumno desarrollará su respuesta así como las explicaciones que considere oportunas. Por favor, razona todas tus soluciones y escribe las explicaciones en azul.
Junto al título de cada apartado se encuentra la puntuación del mismo (pueden obtenerse hasta 10,5 puntos, aunque solamente se evaluará del 0 al 10).
En este apartado, se pide:
library(igraph)
rm(list = ls())
setwd('C:/Users/Javier/Documents/masterAFI/18. Analisis de Grafos/01. Redes sociales/practica/')
dd <- read.csv('red_contactos.csv', sep = ';')
gg <- graph.data.frame(dd, directed = FALSE)
summary(gg)
## IGRAPH ea482fa UN-- 1390 222744 --
## + attr: name (v/c)
vcount(gg)
## [1] 1390
ecount(gg)
## [1] 222744
gg2 <- simplify(gg, remove.multiple = TRUE, remove.loops = TRUE)
# para cada enlace, calculamos cuantos caminos cortos hay entre sus vértices
# el camino mínimo será igual al enlace, por lo tanto, la cantidad será igual al número de enlaces múltiples
# esta cantidad, lo añadimos como peso
E(gg2)$weight = sapply(E(gg2), function(e) {
length(all_shortest_paths(gg, from = ends(gg2, e)[1], to = ends(gg2, e)[2])$res) } )
summary(gg2)
## IGRAPH ea5c390 UNW- 1390 53942 --
## + attr: name (v/c), weight (e/n)
hist(E(gg2)$weight, breaks = 30)
En este apartado, se pide realizar los pasos adecuados para generar un nuevo objeto grafo, que sea conexo, y que involucre a todos los nodos y enlaces de la componente conexa mayor del grafo original.
Comprobamos que no sea conexo:
is.connected(gg2)
## [1] FALSE
Visualizamos los tamaños de las distintas componentes conexas:
ccs <- components(gg2)
ccs$csize
## [1] 1388 1 1
Vemos que hay 3 componentes conexas, una de ellas acumula prácticamente todos los puntos. Nos quedamos con ella:
id_compmayor <- which.max(ccs$csize)
vids <- ccs$membership == id_compmayor
vids <- which(ccs$membership == id_compmayor)
gg3 <- induced_subgraph(gg2, vids = vids )
summary(gg3)
## IGRAPH 1ae7aa5 UNW- 1388 53942 --
## + attr: name (v/c), weight (e/n)
En este apartado, se pide analizar descriptivamente el grafo usando los conceptos que hemos visto durante las clases de teoría:
mean(degree(gg3))
## [1] 77.72622
hist(degree(gg3), breaks = 20)
average.path.length(gg3)
## [1] 3.066029
diameter(gg3)
## [1] 14
Podemos empezar por ver la densidad:
graph.density(gg3)
## [1] 0.0560391
A continuación, ajustamos una Power Law:
deg_dist <- degree_distribution(gg3)
fit <- fit_power_law(deg_dist)
fit
## $continuous
## [1] TRUE
##
## $alpha
## [1] 3.447821
##
## $xmin
## [1] 0.005043228
##
## $logLik
## [1] 391.6632
##
## $KS.stat
## [1] 0.1568648
##
## $KS.p
## [1] 0.03535439
plot(deg_dist + 0.0001, log = 'xy', xlab = 'Node Degree', ylab = 'Probability')
lines(seq(deg_dist), seq(deg_dist)^-fit$alpha, col='red')
Tanto con el p-value como con el plot, podemos ver que nuestros datos no se ajustan a una Power Law.
Podemos tratar de ver que pasaría si eliminamos aquellos nodos con los que obteníamos un valor 0 en la función de distribución:
gg_prueba <- induced_subgraph(gg3, vids = V(gg3)[deg_dist > 0])
deg_dist_prueba <- degree_distribution(gg_prueba)
fit_prueba <- fit_power_law(deg_dist_prueba)
fit_prueba
## $continuous
## [1] TRUE
##
## $alpha
## [1] 3.387724
##
## $xmin
## [1] 0.008536585
##
## $logLik
## [1] 185.4567
##
## $KS.stat
## [1] 0.1363636
##
## $KS.p
## [1] 0.386501
plot(deg_dist_prueba + 0.0001, log = 'xy', xlab = 'Node Degree', ylab = 'Probability')
lines(seq(deg_dist_prueba), seq(deg_dist_prueba)^-fit_prueba$alpha, col='red')
Con el nuevo grafo que hemos creado parece que sí que tendríamos unos datos que se pueden ajustar a una Power Lab.
Podemos observar el número de triángulos por nodo:
triangulos <- count_triangles(gg3)
mean(triangulos)
## [1] 1037.44
hist(triangulos, breaks = 30)
hist(diversity(gg3), breaks = 20, main = '', xlab = 'Diversity')
Vemos las distintas centralidades, hacemos histogramas de ellas y observamos como se correlacionan:
Degree <- degree(gg3)
Eig <- evcent(gg3)$vector
Closeness <- closeness(gg3)
Betweenness <- betweenness(gg3)
par(mfrow=c(2,2))
for (centrality in c('Degree', 'Eig', 'Closeness', 'Betweenness')){
hist(get(centrality), breaks = 30, xlab = centrality, main = '')
}
centralities <- cbind(Degree, Eig, Closeness, Betweenness, triangulos)
round(cor(centralities), 2)
## Degree Eig Closeness Betweenness triangulos
## Degree 1.00 0.93 0.79 0.85 0.92
## Eig 0.93 1.00 0.69 0.78 0.92
## Closeness 0.79 0.69 1.00 0.59 0.59
## Betweenness 0.85 0.78 0.59 1.00 0.87
## triangulos 0.92 0.92 0.59 0.87 1.00
Podemos observar como las centralidades por grado, cercanía y betweeness se parecen bastante entre ellas, además de con el número de triángulos de cada nodo. Sin embargo, parece que Closeness es la que más difiere de las demás.
En este apartado, se pide aplicar dos algoritmos de detección de comunidades, compararlos y seleccionar cuál es, en tu opinión, el que da una mejor respuesta. Razona tu selección.
Usaremos los algoritmos de Infomap y de Walktrap para detectar comunidades, calculamos las modularidades para cada uno y comparamos los resultados de los algoritmos:
comms_infomap <- infomap.community(gg3)
comms_walktrap <- walktrap.community(gg3)
modularity(comms_infomap)
## [1] 0.4686043
modularity(comms_walktrap)
## [1] 0.4634219
compare(comms_infomap, comms_walktrap, method='nmi')
## [1] 0.8131684
Vemos que la función compare nos devuelve un valor de 0.81, es decir, los algoritmos nos devuelven resultados más o menos parecidos. Para elegir uno de ellos, veremos como éstos dividen a sus comunidades:
table(comms_infomap$membership)
##
## 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
## 9 55 392 45 26 11 37 25 12 30 61 87 22 52 61 60 35 4 13 18
## 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40
## 31 35 1 17 38 13 3 25 6 41 10 22 14 3 30 13 7 8 6 6
## 41 42
## 2 2
table(comms_walktrap$membership)
##
## 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20
## 14 12 42 13 134 396 29 50 53 54 28 46 8 15 17 97 2 4 52 83
## 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40
## 5 10 16 7 30 13 19 32 17 2 2 1 1 1 1 1 1 1 1 1
## 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117
## 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
Vemos que el algoritmo de Infomap divide comunidades de una forma más uniforme, mientras que, el algoritmo Walktrap nos devuelve bastantes comunidades formadas por un solo miembro. Con estos resultados, decidimos quedarnos con Infomap como algoritmo de detección de comunidades.
En este apartado, se pide visualizar el grafo coloreando cada nodo en función de la comunidad a la que pertenezca, según tu elección del apartado anterior.
wws <- ifelse(crossing(comms_infomap, gg3), 1, 500)
ll <- layout_with_fr(gg3, weights=wws)
palcol <- rainbow(n=max(comms_infomap$membership))
plot(gg3,
layout = ll,
vertex.label='',
vertex.size=log(degree(gg3)),
vertex.color=palcol[comms_infomap$membership])
Este apartado es el que más peso en la práctica tiene. Vamos a implementar un modelo epidemiológico sobre el grafo que, típicamente, se utiliza para simular escenarios de difusión de enfermedades pero también en contextos como la distribución de rumores e información. Vamos a implementar un modelo SIR que se caracteriza por tener los siguientes parámetros:
Se pide desarrollar una función que tenga como parámetros los tres valores anteriores y un cuarto que sea un grafo que, en nuestro caso, será la componente conexa mayor del grafo original de esta práctica. Dicha función simulará el proceso SIR:
Se pide ejecutar una simulación para tres o cuatro valores del parámetro beta (N y gamma pueden ser fijos en estas simulaciones) de este proceso de manera que se pueda visualizar:
En este primer bloque de código, vamos a definir las siguientes funciones:
generate_sir <- function(N = 5, beta = 0.1, gamma = 0.1, graph){
graph <- set_vertex_attr(graph, 'label', value = 'S')
graph <- set_edge_attr(graph, 'label', value = '')
sample_index <- sample(1:vcount(graph), N)
graph <- set_vertex_attr(graph, 'label', index = sample_index, value = 'I')
n_contagiados = N
n_susceptibles = vcount(graph) - N
n_recuperados = 0
iteraciones_sin_contagios = 0
iteracion= 1
lista_contagiados = c(n_contagiados)
lista_susceptibles = c(n_susceptibles)
lista_recuperados = c(n_recuperados)
lista_nuevos_contagios = c(n_contagiados)
lista_grafo = list(graph)
while (iteraciones_sin_contagios < 3){
n_nuevos_contagios = 0
initial_graph <- graph
for (v in V(initial_graph)){
state <- V(initial_graph)[v]$label
if (state == 'I'){
if (rbinom(n=1, size=1, prob = gamma) == 1){
V(graph)[v]$label <- 'R'
n_recuperados = n_recuperados + 1
n_contagiados = n_contagiados - 1
}
}
if (state == 'S'){
for (neighbor in neighbors(initial_graph, v)){
if ((V(initial_graph)[neighbor]$label == 'I') && (rbinom(n=1, size=1, prob = beta) == 1)){
V(graph)[v]$label <- 'I'
E(graph)[v %--% neighbor]$label <- 'I'
n_contagiados = n_contagiados + 1
n_nuevos_contagios = n_nuevos_contagios + 1
n_susceptibles = n_susceptibles - 1
break
}
}
}
}
if (n_nuevos_contagios == 0){
iteraciones_sin_contagios = iteraciones_sin_contagios + 1
} else {
iteraciones_sin_contagios = 0
}
iteracion = iteracion + 1
lista_contagiados = c(lista_contagiados, n_contagiados)
lista_susceptibles = c(lista_susceptibles, n_susceptibles)
lista_recuperados = c(lista_recuperados, n_recuperados)
lista_nuevos_contagios = c(lista_nuevos_contagios, n_nuevos_contagios)
lista_grafo[[iteracion]] = graph
}
return (list('contagiados' = lista_contagiados,
'susceptibles' = lista_susceptibles,
'recuperados' = lista_recuperados,
'nuevos_contagios' = lista_nuevos_contagios,
'iteraciones' = iteracion,
'grafos' = lista_grafo
))
}
plot_evolution <- function(evolution){
plot(1:evolution$iteraciones, evolution$contagiados,
type = 'l', col = 'red',
ylab = 'Nº Nodos', xlab = 'Iteraciones',
ylim = c(0, max(evolution$susceptibles)),
main = 'Evolución')
lines(1:evolution$iteraciones, evolution$recuperados,
type = 'l', col='green')
lines(1:evolution$iteraciones, evolution$susceptibles,
type = 'l', col='blue')
legend(evolution$iteraciones, max(evolution$susceptibles)/2,
legend=c('Contagiados', 'Recuperados', 'Susceptibles'),
col=c('red', 'green', 'blue'),
lty=1, cex=0.8, xjust = 1)
}
plot_new_infections <- function(list_new_infections, iteraciones){
plot(1:iteraciones, list_new_infections,
type = 'o', col = 'red', log = 'y',
ylab = 'Nº Nuevos Infectados', xlab = 'Iteraciones',
main = 'Nuevos infectados')
}
plot_infection_graph <- function(graph, index){
graph_colors <- ifelse(V(graph)$label == 'I', 'red',
ifelse(V(graph)$label == 'R', 'green', 'blue'))
subg <- subgraph.edges(graph, E(graph)[label == 'I'])
plot(subg,
layout = ll,
vertex.label = '',
vertex.size = log(degree(graph)),
edge.label = '',
vertex.color = graph_colors,
main = paste('Iteración', index))
}
create_gif <- function(list_graphs){
img <- magick::image_graph(width = 1000, height = 1000)
for (i in 1:length(list_graphs)){
plot_infection_graph(list_graphs[[i]], i)
}
dev.off()
animation <- magick::image_animate(img, fps = 2, optimize = TRUE)
print(animation)
}
Veremos distintas ejecuciones del simulador del proceso SIR, dibujando las gráficas y los GIF correspondientes a cada una de ellas.
Nuestras distintas ejecuciones tendrán como parámetros los mismos valores para gamma (probabilidad de recuperación; gamma = 0.1), N (número de nodos iniciales infectados; N = 5) y el grafo (utilizaremos el grafo resultante de los apartados anteriores). El parámetro que vamos a ir variando será el beta.
ejemplo_SIR_1 <- generate_sir(graph = gg3, beta = 0.1)
plot_evolution(ejemplo_SIR_1)
plot_new_infections(ejemplo_SIR_1$nuevos_contagios, ejemplo_SIR_1$iteraciones)
create_gif(ejemplo_SIR_1$grafos)
## format width height colorspace matte filesize density
## 1 gif 1000 1000 sRGB TRUE 0 72x72
## 2 gif 784 859 sRGB TRUE 0 72x72
## 3 gif 442 687 sRGB TRUE 0 72x72
## 4 gif 605 809 sRGB TRUE 0 72x72
## 5 gif 620 857 sRGB TRUE 0 72x72
## 6 gif 616 792 sRGB TRUE 0 72x72
## 7 gif 620 857 sRGB TRUE 0 72x72
## 8 gif 613 857 sRGB TRUE 0 72x72
## 9 gif 627 857 sRGB TRUE 0 72x72
## 10 gif 592 857 sRGB TRUE 0 72x72
## 11 gif 590 857 sRGB TRUE 0 72x72
## 12 gif 535 856 sRGB TRUE 0 72x72
## 13 gif 587 859 sRGB TRUE 0 72x72
## 14 gif 546 857 sRGB TRUE 0 72x72
## 15 gif 769 856 sRGB TRUE 0 72x72
## 16 gif 595 687 sRGB TRUE 0 72x72
## 17 gif 754 695 sRGB TRUE 0 72x72
## 18 gif 769 856 sRGB TRUE 0 72x72
## 19 gif 244 682 sRGB TRUE 0 72x72
## 20 gif 266 738 sRGB TRUE 0 72x72
## 21 gif 490 653 sRGB TRUE 0 72x72
## 22 gif 349 688 sRGB TRUE 0 72x72
## 23 gif 769 857 sRGB TRUE 0 72x72
## 24 gif 241 687 sRGB TRUE 0 72x72
## 25 gif 750 781 sRGB TRUE 0 72x72
## 26 gif 206 630 sRGB TRUE 0 72x72
## 27 gif 426 578 sRGB TRUE 0 72x72
## 28 gif 205 685 sRGB TRUE 0 72x72
ejemplo_SIR_2 <- generate_sir(graph = gg3, beta = 0.05)
plot_evolution(ejemplo_SIR_2)
plot_new_infections(ejemplo_SIR_2$nuevos_contagios, ejemplo_SIR_2$iteraciones)
create_gif(ejemplo_SIR_2$grafos)
## format width height colorspace matte filesize density
## 1 gif 1000 1000 sRGB TRUE 0 72x72
## 2 gif 784 859 sRGB TRUE 0 72x72
## 3 gif 443 786 sRGB TRUE 0 72x72
## 4 gif 612 809 sRGB TRUE 0 72x72
## 5 gif 626 792 sRGB TRUE 0 72x72
## 6 gif 613 792 sRGB TRUE 0 72x72
## 7 gif 619 809 sRGB TRUE 0 72x72
## 8 gif 627 792 sRGB TRUE 0 72x72
## 9 gif 614 856 sRGB TRUE 0 72x72
## 10 gif 620 856 sRGB TRUE 0 72x72
## 11 gif 614 857 sRGB TRUE 0 72x72
## 12 gif 613 856 sRGB TRUE 0 72x72
## 13 gif 615 856 sRGB TRUE 0 72x72
## 14 gif 599 856 sRGB TRUE 0 72x72
## 15 gif 592 857 sRGB TRUE 0 72x72
## 16 gif 592 856 sRGB TRUE 0 72x72
## 17 gif 562 695 sRGB TRUE 0 72x72
## 18 gif 554 856 sRGB TRUE 0 72x72
## 19 gif 574 702 sRGB TRUE 0 72x72
## 20 gif 592 857 sRGB TRUE 0 72x72
## 21 gif 592 856 sRGB TRUE 0 72x72
## 22 gif 594 856 sRGB TRUE 0 72x72
## 23 gif 592 856 sRGB TRUE 0 72x72
## 24 gif 286 681 sRGB TRUE 0 72x72
## 25 gif 343 639 sRGB TRUE 0 72x72
## 26 gif 592 856 sRGB TRUE 0 72x72
## 27 gif 614 857 sRGB TRUE 0 72x72
## 28 gif 524 854 sRGB TRUE 0 72x72
## 29 gif 615 857 sRGB TRUE 0 72x72
## 30 gif 467 694 sRGB TRUE 0 72x72
## 31 gif 578 738 sRGB TRUE 0 72x72
## 32 gif 361 692 sRGB TRUE 0 72x72
ejemplo_SIR_3 <- generate_sir(graph = gg3, beta = 0.01)
plot_evolution(ejemplo_SIR_3)
plot_new_infections(ejemplo_SIR_3$nuevos_contagios, ejemplo_SIR_3$iteraciones)
create_gif(ejemplo_SIR_3$grafos)
## format width height colorspace matte filesize density
## 1 gif 1000 1000 sRGB TRUE 0 72x72
## 2 gif 784 859 sRGB TRUE 0 72x72
## 3 gif 441 628 sRGB TRUE 0 72x72
## 4 gif 442 676 sRGB TRUE 0 72x72
## 5 gif 473 693 sRGB TRUE 0 72x72
## 6 gif 565 691 sRGB TRUE 0 72x72
## 7 gif 545 693 sRGB TRUE 0 72x72
## 8 gif 606 809 sRGB TRUE 0 72x72
## 9 gif 600 792 sRGB TRUE 0 72x72
## 10 gif 619 792 sRGB TRUE 0 72x72
## 11 gif 620 792 sRGB TRUE 0 72x72
## 12 gif 617 792 sRGB TRUE 0 72x72
## 13 gif 620 792 sRGB TRUE 0 72x72
## 14 gif 603 793 sRGB TRUE 0 72x72
## 15 gif 613 792 sRGB TRUE 0 72x72
## 16 gif 612 792 sRGB TRUE 0 72x72
## 17 gif 620 792 sRGB TRUE 0 72x72
## 18 gif 612 792 sRGB TRUE 0 72x72
## 19 gif 613 792 sRGB TRUE 0 72x72
## 20 gif 604 792 sRGB TRUE 0 72x72
## 21 gif 587 809 sRGB TRUE 0 72x72
## 22 gif 603 792 sRGB TRUE 0 72x72
## 23 gif 620 792 sRGB TRUE 0 72x72
## 24 gif 539 792 sRGB TRUE 0 72x72
## 25 gif 603 792 sRGB TRUE 0 72x72
## 26 gif 603 792 sRGB TRUE 0 72x72
## 27 gif 568 792 sRGB TRUE 0 72x72
## 28 gif 603 792 sRGB TRUE 0 72x72
## 29 gif 520 710 sRGB TRUE 0 72x72
## 30 gif 603 792 sRGB TRUE 0 72x72
## 31 gif 532 703 sRGB TRUE 0 72x72
## 32 gif 560 703 sRGB TRUE 0 72x72
## 33 gif 604 792 sRGB TRUE 0 72x72
## 34 gif 560 785 sRGB TRUE 0 72x72
## 35 gif 435 652 sRGB TRUE 0 72x72
## 36 gif 512 681 sRGB TRUE 0 72x72
## 37 gif 393 682 sRGB TRUE 0 72x72
## 38 gif 271 602 sRGB TRUE 0 72x72
## 39 gif 384 558 sRGB TRUE 0 72x72